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O . Abstract: An algorithm is presented here, for discovering Hopf-Bifurcation varieties of 

polynomial dynamical systems. It is based on the expression of specific polynomials, as 
sums of products of first degree polynomials, with parametrical coefficients. By giving to 
these parameters certain values, we ensure the positiveness of some quantities, construct- 
ing thereby proper Lyapunov functions, which guarantee the stability of the equilibrium 
point. The points where the afore mentioned positiveness fails, define the Hopf-Bifurcation 

/-N , varieties upon discussion. 
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1 Introduction 

It is widely known that when engineers and economists have to analyze mechanical, electri- 
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O . cal or economical dynamic phenomena, they are usually dealing with nonlinear dynamical 






systems. They give attention to the stability of those systems and to some special dynamic 

properties they possess (the existence of limit cycles, for instance), p. Another issue of 

O ' practical importance is whether the systems maintain their dynamic behavior as certain 

^ ■ parameters are varied. Especially, it is interesting if a given equilibrium point remains 

stable or unstable. These parameters are called bifurcation or Hopf-Bifurcation parameters 

and the values at which changes occur, if any, are called Hopf-Bifurcation points, [5].|1U]. 

r> I The aim of this paper is to present an algorithm which detects Hopf-Bifurcation points of 

c^ ' autonomous polynomial dynamical systems of the form: 

i=$(x,^) (1) 

where, x is the state vector, consisting from functions of t, $ is a vector function, consisting 
from polynomial functions of elements of x and /i a set of parameters. Actually, the said 
algorithm, discovers the so called bifurcation varieties or bifurcation curves of a given 
equilibrium point of ([T]). These are sets of values of the parameters which satisfy certain 
relations. The violence of those relations implies the change of the stability behavior of the 
system around the equilibrium point. 

There are a lot of efforts in the literature, for the description of proper algorithms which 
help us to determine bifurcation points. Let me refer to [6] , [7j , [TT] , to mention but a few. 
The algorithm presented here works as follows. First it accepts a polynomial function L 
as a Lyapunov function candidate, for the system ([T]), and an equilibrium point xq. Then 
it calculates the derivative of L across the trajectories of ([T]), denoted by L. After that, it 



checks the positiveness of the quantities L and —L. If this is true, under the assumption 
that the parameters satisfy a certain set of relations denoted by B, then we have stabihty of 
the equihbrium point and the boundary of B defines the Hopf-Bifurcation variety. Indeed, 
" crossing " this boundary it means that the relations B do not hold any more and the 
stability of xq collapses. 
To explore the positiveness of the quantities L and — L, we decompose them as follows: 

V = Ci(M^,,.,^)[iyi,_i,i+Xi]^'^-^-[M^2,-l,l + W^2A,l^l+^2r^-[W^3,-l,l + W^3,l,i:^^ ■ 

■ ■ ■ [M^n,-l,l + Wn,l,lXi + 1^,2,1X2 + " " " + X„]^"'^ + 
+C2(l^i,a,^)[m,-l,2 + Xlp'^-[H^2,-l,2 + H^2,l,2Xl+X2]^'^-^-[K^3,-l,2 + H^3,l,2Xl + V^^^^^ 

■ ■ ■ [Wn,-1,2 + W;,l,22:i + Wn,2aX2 + " " " + Xn]^"'^ + ■ ■ ■ 
+Cfc(W.,.,^)[M^l,_l,fc+Xi]^'i-*-[H^2,-l,fc + H^2,l,fcXl+X2]^'^-*-[H^3,-l,fc^ 

■ ■ ■ [Wn,-l,k + Wn,l,kXl + Wna,kX2 + " " " + Xn]^"'" + Rw (2) 

where the exponents ja,b are specific positive whole numbers. The quantities Wi^a,^ are un- 
determined parameters that can take real values, the coefficients Cj(VFi,o-,(/.) are depending 
on the parameters VFj,o-,0 and the quantity Ry^ is a polynomial of the parameters Wi^^j^^i, 
only, called the remainder. We obtain this " factorization " of the polynomials by means of 
a recursive algorithm, introduced in [2], which resembles to the Euclidean Algorithm, and 
annihilates successively the maximum terms. Then, we seek for those values of the param- 
eters Wi^cr,<j>, which eliminate the non-square terms and make the coefficients of the square 
terms and the remainder, positive. Obviously, if this can be achieved, the positiveness of L 
and —L is secured. 
The main advantages of the method are: 

1) Our approach is symbolic in nature and not numeric. 

2) The calculations can be easily carried out, since the coefficients in the expression 
([2]) have a specific construction. Each of them contains a number of parameters which 
is larger or equal than the number of the parameters of the previous coefficient. This 
triangular structure, more known as a sparse system of algebraic equations, permits their 
easier handling, |12j . 

3) The bifurcation values of the parameters fi, can be found straightforward, since they 
are involved in the calculations together with the artificial parameters Wj,o-,0, and therefore 
can be considered as polynomial functions of them. 

We have to make clear that our method does not provide necessary and sufficient conditions. 
In other words, if our approach fails this does not mean that there are not bifurcation 
varieties or that another method could not find them. 
Throughout this paper R will denote the set of real numbers. 

2 Preliminaries 

In this section we present the basic tools on which the basic algorithm is relied. Let 
(xi,a;2, . . . ,Xn) be a vector of rz-indeterminates which abbreviate by x. An expression of 
the form p = J2x=i c\Xi^'^X2^'^ ■ ■ ■ Xn"'^ , where Ca G R and some of the exponents aij E Z"*" 
are not equal to zero, is called a polynomial in 2;i,a;2, . . . ,Xn with real coefficients or, for 



short, a real polynomial. An element x^'^x"^'^ ■ ■ ■ Xn"'^ is called a monomial and an element 
cxXi^'^X2^'^ ■ ■ ■ Xn"'^ is called a term. The quantity cx is the coefficient of the term. The sum 
fli,A + fl2,A + ■ ■ ■ + «n,A is Called degree of the term. A term is called even if all of its exponents 
are even integers, otherwise it is called an odd term. The term which corresponds to the 
exponent (0, 0, . . . , 0), is the constant term. If we use the multi-index notation slx G (Z"'")"' 
to denote the vector a^ = {aix, 0,2^X1 ■ ■ ■ , fln,A), we write a monomial compactly as x*^^ and 
a polynomial as p = J2'x=i cax^^. The set of all real polynomials in xi, ^2, . . . , a;„ is denoted 
by R[xi, X2, . . . , Xn] or R[x]. Let (pn,\ = Xi^'^X2^'^ ■ ■ ■ Xn"'^ and 0m,/^ = Xi^'''x2^''' ■ ■ ■ x^-*" 
be two monomials. We define the lexicographical order among monomials [1], as follows: 
we say that (j)n,\ is ordered less than (pm,^, denoted by 0n,A -< (pm,fi, if either n < m or 
n = m and in the vector difference 0„^^ — 0m,A the left-most nonzero entry is positive. 
In other words, the monomials are ordered as follows: Xi -<■■■-< x[ -<■■■ ^ X1X2 -< 
■ ■ ■ -< X1X2 -< ■ ■ ■ ^ X1X2X3 ^ ■ ■ • Let p be a given polynomial, ordered lexicographically, 
the term that corresponds to the maximum monomial is called the maximum term denoted 
by maxterm{p), its degree is called the polynomial degree and it is denoted by deg(p, x). 
The next definitions will play a crucial role in the subsequents. 

Definition 2.1 Let 7rj(x), i = 1, . . . ,m be a collection of polynomials in R[x]. Then we 
set 

V = V{TVi) = {9eK' : 7ii{e) = 0, for all i = 1, . . . ,m} 

We call V the variety defined by 7rj(x), i = 1, . . . ,m. 

Definition 2.2 Let 7rj(x), i = 1, . . . ,m be a collection of polynomials in R[x]. Then we 
set 

A = A(7r,) = {6' G R" : TTi{e) < 0, for all i = 1, . . . ,m} 

We call A the semi-algebraic set defined by 7ri(x), i = 1, . . . , m. 

Definition 2.3 Let A be a semi-algebraic set defined by the polynomials 7rj(x), i = 1, . . . ,m. 
The variety {6 G R" : vrj(0) = 0} is called the boundary of A and it is denoted by dA. 

Obviously, if we reverse the sense of the inequalities, nothing will change in the meaning of 
the above definitions. 

Let us have a polynomial p = J2x=i CA(ct)x**^; where the coefficients are polynomial expres- 
sions of a set of certain parameters q = {qi,q2, ■ ■ ■ ,qk)- Let us further suppose that we have 
the variety: 

^7 = {r = (ri,r2, . . . ,rfc) G R' : 7r,(r) = 0, ^ = 1, . . . , ^} C R*= 

where ttj is a collection of polynomials of q. We say that the polynomial p is evaluated over 
the variety U, thus writing p\u, if the parameters q take values from the set U. Rigorously, 
we have: 

^It/ = I E CA(r)x^\ reu\ 

If ?7 is a finite set then p\u is a finite set too, consisting from polynomials of x, with real 
coefficients. If U is an infinite variety, then p\u is an infinite set. If, furthermore, in this 
case a parametrization for the description of the variety is available, we can use it for the 



description of the set p\u, too, in a natural way. Indeed, let us assume that we have the 
real functions (pk '■ R'*' — >■ R, A < .^ — 1 and that the points given by the relations 

rk = iPk{to,ti-,---,t\), A<^-1, A; = 0, ...,^- l,to,...,tA G R 

lie in U . These functions constitute a parameterization for the variety U and the parameters 
q. Then, 

[a=o 
The above terminology can be extended to a set C of polynomials, as follows: 

^u = U P\u 

Let us now suppose that we have n functions of t : xi(t), X2(t), . . ., x„(t) and n polynomial 
functions of Xi,X2,...,x„, i.e., <l>i(xi, X2, . . . , x„), $2(a;i, X2, . . . , x^), . . ., $„(xi, X2, . . . , x„). 
An expression of the form: 

Xi = ^i{xi,X2,...,Xn) 

X2 = ^2{,Xi,X2,...,Xn) (3) 



Xn ^n\Xij X2i ■ ■ ■ i Xji) 

where by Xk we denote the derivative of the Xk function, with respect to the time, it is called 
an autonomous polynomial dynamical system. We write it compactly as x = $(x) where 
x(t) = {xi(t), X2{t), . . . , x„(t)) G R", for each t, is the state space vector. If the coefficients 
of the polynomials $j depend on a certain set of parameters fi = (/ii, /U2, . . . , Hk) we write 

x= $(x,^) 

A state xq is called equilibrium point of (j4]), if and only if $(xo) = or $(xo, /i) = 0. In 
the latter case xq may be depend on the parameters /i. 

The next notions are classical in the literature, but we include them here for the sake of 
the self-reliance of the paper, [5j. 

Definition 2.4 An equilibrium point Xq of ^ is stable if for each e > there is S = S{e) > 
such that 

||x(0)|| <(5^ ||x(t)|| <e, Vt>0 

Definition 2.5 The equilibrium point xq is unstable if it is not stable 

Definition 2.6 The equilibrium point xq is asymptotically stable if it is stable and 6 can 
be chosen such that 

||x(0)||<(5^ hm ||x(t)||=0 

t — ^ + 00 



To determine stabilizability of a given equilibrium point, we adopt the Lyapunov function 
methodology. Let L : D —>■ H, he a. continuously differentiable function defined in a domain 
D C R" that contains the equilibrium point. The derivative of L along the trajectories of 
([3]), denoted by ^(x), is given by: 

We are ready now to establish the Lyapunov's stability theorem. 

Theorem 2.1 f^ Let xq be an equilibrium point for ^ and D C R" be a domain con- 
taining xq. Let L : D ^ Ti be a continuously differentiable function such that 

L(0) = and L(x) > m D- {xq} 

L(x) <0 m D 
Then xq is stable. Furthermore, if 

L(x) <0 %n D 

then Xq is asymptotically stable. 

The next definition introduces the notion of Hopf-Bifurcation variety (or Hopf-Bifurcation 
curve), which is a generalization of the well known Hopf bifurcation point. It is nothing else 
than a set of relations for the parameters /i the violation of which can make the equilibrium 
point either stable or unstable, depending from the kind of violation. Specifically: 

Definition 2.7 Let x = $(x, /i) be an autonomous polynomial dynamical system depending 
from a set of parameters /i = (/ii, H2, ■ ■ ■ , fJ'k) o-nd xq be an equilibrium point. Let, moreover, 
7rj(/i), J = 1, . . . , z/, be a collection of polynomials of fi. Let us have the semi- algebraic set 
S = {9 E R'^ : TijiO) < 0}. //xo is asymptotical stable, for the system x = $(x, /z)!^ then 
the set B = dS is called the Hopf-Bifurcation variety of the parameters /i, at the point xq. 



Example 2.1 Let us have the system 

X = —yx + kx 



y = ux^ 

We write it compactly as yi = $(x, /x), with x = {x,y) and fj, = {k^p). This system has 
an infinite number of equilibrium points of the form (0, A), A G R. If S and M. are the 
varieties S = {k = l^v = 1} , M. = {(A;, ly) : k = t,iy = 2t,t E R}, then the evaluation of 
the system over these varieties will give 

i=$(x,/.)i^^ •17"^^"' 

y '*' 
x = $(x,/.)U^ ^ = 2tx2 



Let us now check if L = vx^ + {y — A)^ is a Lyapunov function. Indeed, L(0, A) = and 
L(x) > for u > 0. Furthermore, i/(x) = 2i/{k — A)x^, which is negative or equal to zero 
if k < X. Therefore, the evaluation of the original system over the semi- algebraic set 

S = {{k,u) : A; < A,z/ > 0,A G R} 

will provide us with a dynamical system, which is either stable or asymptotical stable (de- 
pending on the fact being k = X or k < X) at the equilibrium point (0, A). The boundary of 
S is dS = {(A;, z/) : A; = A, z/ = 0, A G R} and it consists a bifurcation variety. 

3 The Algorithm 

In this section we present the synibohc algorithm which discovers bifurcation varieties. This 
algorithm annihilates step by step the current maximum terms, by subtracting a suitable 
product of first degree polynomials with parametrical coefficients. We suppose that two 
algorithms are available, the first, named Solve- Algorithm, solves a system of polynomial 
equations {p = : p G P}, the second, named InSolve- Algorithm, solves a system of 
inequalities {p < : p G P}. The construction of such algorithms are the subject of the 
research of computational algebra, and certain methodologies have been developed toward 
this scope, [1],[3],[1]. 

THE FQRMAL-BIF-ALGORITHM 

Input: The integer m, the polynomials <I>j(x, /x), i = 1, . . . , /, the point xq = (xo,i, xo,2, • • • , xo^n), 
the sets of parameters W = {I^i,^,^}, S = {Sij^k} and A = {Aj-^,J2, .-■,«*} 

Output: The set B. 

Step 1: Construct the polynomial ( Lyapunov Function Candidate). 

m 
^ — Z^ Z^ ^(«l,«2v,ifc)^l ^2 ■ ■ "^fc 

fc=i(Ji,i2,...,ife)e/fcCZ+* 

where ^(j^^ij, •••,«*:) ^ A. If 1^ = , then the corresponding terms and their coefficients 
do not exist. 

Step 2: Use the Formal-Subroutine[L, S], that is with inputs the polynomial L and the 
parameters S, to obtain the sets Ol,El. The last element of El is denoted by Rl- 



Step 3: Construct the polynomial 



i=l 



dxi 



Step 4: Use the Formal-Subroutine[y, W], that is with inputs the polynomial V and the 
parameters W, to obtain the sets OvjEy. 

Step 5: Set O = O^ U Oy U {L{xo)}, E = ElU Ey, 



Step 6: By means of the Solve-Algorithm, find the values of the parameters A, fi, S, 
W for which G = 0, for all G G O . They construct a variety denoted by O. 

Step 7: Create the semi-algebraic set: 

J = {(A, ^, S, W) G R'^ : G > 0, G G {El - {Rl})\o, Rl\o >0,H > 0,H g Ev\o} 

Step 8: By means of the inSolve-Algorithm check the feasibility of J, that is if J / 
or not. 

Step 9: IF J = THEN the method fails, stop ELSE denote by K the set of 
equations which defines the boundary dJ and then, set 

B = K-{eeK, with, deg(Q,/i) = 0} 

Goto the output. 

THE FORMAL SUBROUTINE [p,W] 

Input: A set of undetermined parameters W = {Wi^a,ip}, taking values in R. 
A multivariable polynomial 

e 
p=J2cx{f^)xT''---x>' , ai,AGNU{0}, i = l,...,n, X=l,...,e 
A=l 

where the coefficients cx{fi) depend from the vector of parameters /i = (^i, . . . ,^,y). 

Output: The sets Ep and Op. 

Initial Conditions: k = 0, Rq = p, Ep = {}, Op = {} 

REPEAT UNTIL Rk does not depend on any of the variables xi,X2,..., x„, THEN 

set Ep = EpU{Rk}. 

Step 1: Setk = k + 1. 

Step 2: Find the maximum term of Rk-i, maxterm{Rk-i) = Cfc(//, W)x'{^''' •••x^"''°. 
The coefficient Cfc(/i, W), in the first iteration, is either a constant number or de- 
pends from the parameters /x only. Then, it depends on the set of parameters W, 
too. 

Step 3: IF at least one of the exponents jjfe, z = 1, . . . ,n is an odd positive integer 
THEN Op = Op\J {cfe(/i, W)} ELSE Ep = EpU {cfe(;U, W)} 

Step 4: Form the polynomials: 

L2,k = W2,^l,k + W2,l,kXl + X2 



Ln,k = Wn-l,k + Wn,l,kXl + Wn,2,kX2 + Wn,Z,kX^ H V Xr. 



Step 5: Make the subtraction: 

RETURN 

The next two theorems describe the behavior of the Formal-Bif-Algorithm. 

Theorem 3.1 The Formal-Bif-Algorithm terminates, after a finite number of steps. 

Proof: The termination of the algorithm depends from the termination of the Repeat 
procedure which appears at the Formal- Subroutine. But it can be easily proved that the 
operation at step 5 of the subroutine, annihilates the current maximum term of Rk-i and 
thus, the maximum term of Rk will be ordered less accordingly to the lexicographical order. 
Therefore, eventually, all the terms which contain at least one of the variables Xi,X2, ■ ■ ■ ,Xn 
will be eliminated and the procedure will stop. 

Theorem 3.2 The amount of the undetermine parameters Wijk, appeared at the coefficient 
Cjt(/i, W) of the step 5 of the subroutine, is larger or equal to the amount of parameters Wijk, 
appeared at the coefficient Cfc_i(/i, W). 

Proof: Let px^x^^ ■ ■ ■ x^ be the maximum term of the polynomial p. When we visit step 
5 of the subroutine for first time, we get Ci(/i, W) = p. This coefficient contains a zero 
amount of parameters Wijk- The next higher ordered term is Aa;;^^~ X2^ ■ ■ ■ x'^"- This term 
may exist at the polynomial p but it will be also created by the product L^^LaJi ■ ■ ■ L^"i= 
(VFi^_i^i + xi)'^'^{W2-i^k + W^2,i,fc2^i + 3^2)^^" ■ ■ of the step 4 of the subroutine. Using the 
binomial theorem (Newton's theorem) for the expression (W^i _i 1 + Xi)^'^, we finally get: 
C2(/U,W) = A — hiWi^-.i,i. This coefficient contains obviously a larger amount of W- 
parameters than the previous one. For the coefficient C3(^, W) we have C3(/i, W) = g — 
h2Wi-ifl — ^^ '^~ Wl_i I, which has the same property and therefore, working inductively, 
we can establish the theorem. 

The next theorem is devoted to the efficiency of the algorithm to the discovering of the 
bifurcation varieties. 

Theorem 3.3 Let us suppose that we have the nonlinear dynamical system 

i = $(x,^) (4) 

where $ = ($1, $2, • • • , "^p) o-i^^ polynomials of ^, with polynomials expressions of fi, as 
coefficients, and xq an equilibrium point. We apply the Formal-Bif-Algorithm and let B be 
its output. If B ^ ^, then all the members of B define the H op f- Bifurcation variety of Q), 

at Xq. 

Proof: We shall follow the Formal-Bif-Algorithm step by step. Step 1 will create a Lya- 
punov function candidate. Step 3 will create the quantity V , which is the opposite of its 
derivative across the orbits of the system. Both of them are polynomials of x with para- 
metrical coefficients. Applying the Formal-Subroutine we get the next expressions for L 
and V correspondingly: 

k 
T — \^ ? T^^-'' T^^'^ . . . 7'J'".'- -I- P. 



r=l 



V — T, CrL^'j! L2r ■ ■ ■ L^j!' + iiv 



with 



and 



-^2,r = '-'2,-l,r + ^2,l,rXl + X2 

^l,r = W^l,-l,r + a^l 

^2,r = W2-l,r + ^^2,1,^3:^1 + ^2 



The coefficients 5^,0^ are polynomial expressions of the parameters A,W, S and /i. We 
obtain these expressions by backward substitution of the values of i?^, obtained in the step 
5 of the subroutine. Step 6, in combination with the construction of the set O, finds those 
values of the parameters A, //, S and W, which guarantee that all the odd terms of L and V 
will be eliminated. Furthermore they ensure that the Lyapunov function candidate is equal 
to zero at the equilibrium point xq. Step 7 of the algorithm, constructs a set of variables 
which make the coefficients of the even terms of the L strictly positive, the remainder R^ 
positive or equal to zero and the coefficients of the even terms of V strictly positive or zero. 
All the above indicate that the above expression L is a Lyapunov function for the nonlinear 
system, which is negative across its orbits and thus xq is stable or asymptotical stable 
depending if y < or V < 0. Step 9 provides us with the boundary of the previous set, in 
other words this particular set which "separates" the two different " behaviors " those of 
stability and instability of the equilibrium point. The subset of the above set consisting only 
with equations which involves the parameters /i, (step 9), defines the bifurcation variety. 

Remark 3.1 The keynote of the whole process is the solvability of the polynomial equations 
or inequalities, appeared at the steps 6 and 7 of the algorithm. This is the classical problem 
in algebraic geometry and certain methods have been developed toward this direction. It still 
remains a hard problem. Nevertheless, the polynomials appeared at the sets O or J have a 
particular structure, which make the solution of the equations easier. Indeed, accordingly to 
theorem lS.Si they have a "triangular" construction and can be faced via methods of solution 
devoted to sparse systems, JJl 



4 Examples 

To clarify the previous algorithm we present certain examples. 

Example 1. To indicate how the algorithm works in practice, we deal ffistly with a rather 

simple example. We consider the system 

LljtAJ iXj 



y = -y 

In this case x = (x, y), // = /i and $ = ($i, $2) with $1 = /ux — x^ and $2 = —y- This 
system has two equihbrium points (0,0) and (/i,0). By means of classical tools, that is the 
Jacobian, we can show that (0,0) is a stable point ioi fi < and a saddle point for /i > 0, 
and that (/i, 0) is a stable point for /i > and a saddle point for /i < 0, [5J. 
Let us handle this system by using the algorithm developed previously. We shall work 
with the point (0,0). We take as L a second degree homogeneous polynomial, that is 
L = Aix^ + ^2?/^ + ^sxy, and we follow the Formal-Bif-Algorithm step by step. By 
substituting backwards the results of the step 5 of the subroutine we shall take the next 
expression for the polynomial L: 

L = A2{S2,-i,i + S2,i,ix + yf + (A3 - 2^2^2,1,1) (^1,-1,2 + a;)(52,-i,2 + ^2,1,23; + y)+ 
+(-y43S'i _i,2 - 2^25*2-1,1 + 2^25*1 _i,25'2,i,i) (5*2 -1,3 + S'2,i,3a; + y)+ 

+ (Ai - ^2^2' 1,1 - ^3>^2,l,2 + 2A252,l,l'52,l,2)(5l,-l,4 + xf + ■ ■ ■ + Rl 

with 

Rl = ~^i*S'i,_i,4 + 2yiiS'i,_i_4S'i,_i_5 — A2S2_ii + ■ • • 

For the quantity V = 2Aix^ — 2^Aix'^ + 2A2y'^+A^xy + A^x^y — fiA^xy, we get: 

V = 2^2(1^2,-1,1 + W2,i,ix + yf + ^3(1^1,-1,2 + xf{W2,^i,2 + W2,i,2X + y)+ 
+ (-31^1,-1,2^3) (Vri,-i,3 + a^)'(1^2,-i,3 + 1^2,1,32^ + y) + 

+ (A3-/iA3-3A3ll^i%i,2 + 6A3Vl'l,-l,2Vl^l,-l,3-4A2ll"2,l,l)(Vl^l,-l,4 + 2:)(iy2,-l,4 + 11^2,l,4a; + ?/) + 

+ci(1^2,-i,5 + 11^2,i,5X + y) + (2Ai - 1^2,1,2^3) (W'l.-i.e + xY+ 

+C2(m,-1,7 + xf + C3(l^l,-l,8 + xf + C4(l^l,-l,9 + x) + Ry 

(We do not write the coefficients ci, C2, C3, C4 and the remainder i?y, exphcitly due to their 
large size). Thus, the set O, consisting from all the coefficients of the odd terms of L and 
y, is 

= {v43- 2^2^2,1,1, -^351,-1,2- 2^2^2,-1,1 + 2v425l,-l,252,l,l,..., ^, -31^1,-1,2-43, 
A3 - M3 - 3A3l^i%i,2 + 6^31^1,-1,21^1,-1,3 - 4^21^2,1,1, • • • , 0} 

The last number 0, corresponds to the evaluation of the Lyapunov function candidate at 
the equilibrium point xq = (0,0). The sets El and Ey, consisting from the coefficients of 
the even terms of L and 1/, correspondingly, are: 

El = {^2, Ai — ^25*2^1 1 — ^35*2, 1,2 + 2v42S'2, 1,1 5*2,1,2, ■ ■ ■ , Rl} 

Ev = {2A2, 2Ai - W2,i,2A3, . . . , Rv} 

Now, the variety O, consisting from the values of the parameters which vanish the non- 
square terms, is: 

= {^2,-1,1 = 0, 11^2,-1,1 = 0, Si,_i,4 = 0, 1^2,1,1=0, ^2,1,1=0, 

10 



Evaluating El and Ey over the values of O, we get: 

El\o = {AuA2} , Ev\o = {2Ai,-2Aifi,2A2} 

These values of the coefficients produce the next expressions for L and V: 

L\o = Aix^ + A2?/^ V\o = 2Aix^ - 2Aiiix^ + 2^2?/^ 

Now, the semi-algebraic set J is 

J = {(Ai, A2, ^),Ai> 0, A2 > 0, 2Ai > 0, -2Ai^ > 0, 2A2 > 0} 

which is feasible for v4i > 0, j42 > 0, /i < 0. In other words, the above values guarantee 
the positiveness of L and V, which ensure the stability of the origin. The boundary dJ 
is defined by the equations Ai = 0,A2 = 0, fi = 0. The only equation which involves the 
parameter is /i = 0, this equation defines the bifurcation variety B = {/i : /j, = 0}. This 
result coincides with that provided by the classical theory, [5]. We can repeat a similar 
analysis for the point (/i, 0), too. 

Example 2. In this example we exhibit the applicability of the method in the case of one 
dimension dynamical systems. Let us consider the system: 

X = ax^ + bx'^ + ex + d 

To simplify the manipulation we take d = —{a + b + c). This means that x = 1 is an 
equilibrium point. We shall work with this specific point. As Lyapunov function candidate 
we shall use the quantity L = x"^ + Aix + 1. The Formal Bif- Algorithm will give: 

L={x + Si,_i,i)2 + {A, - 25i,_i,i)(x + 5i,_i,2) + 

and 

V = -2a{x + 1^1,-1,1)^ + (-26 - aAi + 8aW^i,_i,i)(x + ^,-1,2)'+ 

+Ci(x + VI/1,-1,3)^ + C2{X + 1^1,-1,4) + Rv 

We seek those values of the parameters, which will eliminate the non-even terms. This is 
achievable ii b = —3a, and thus O = {S'l _i_i = —1, Ai = —2, PFi _i_3 = —1, VFi^_i^i = —1}. 
Evaluating L and V over O we take: 

L\o = ix-lf , V\o = i6a-2c){x-lY-2a{x-iy 

The positiveness of V is guaranteed by the feasibility of the set J = {(a, c) : 6a — 2c < 
0, — 2a > 0}. This set is indeed non-void and dJ is defined by the equations 69 — 2c = 0, 
for any ^ < and a = 0. Thus, the bifurcation variety is {{a,c) : a = 0,c = 3^},for any 
given 6 < 0. 

Example 3. To illustrate the compatibility of the method with linear systems, we examine 
the next case: 

X = ax + /3y 
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The origin (0, 0) is the only equihbrium point. Let us work with the function L = x^ + -ip. 
This is a classical Lyapunov function and thus we shall not deal with it any more. For the 
quantity V the Formal Bif- Algorithm will give: 



V 



= -2S{W2,-is + W2,i,ix + yf + (-2/3 - 27 + 45l^2,i,i)(W^i,-i,2 + x) (1^2,-1,2+ 

+W2,i,2X + y) + m,^i,2(2/9 + 27 + 45 - A5W2,i,i) 
(W^2,-i,3 + W2,i,^x + y) + [-2a + 25W2,i,i{W2,i,i - 2^^2,1,2) + 2W^2,i,2(/? + 7)] 



■(W^l,-l,4 + xf + Ci(W^i,„i,5 + x) + /2v 



The next set of values of the parameters O = {<S'2,-i,i 
W^i,-i,4 = 0, ^2,1,1 = 0, 1^2,1,1 = ^}, Will give: 



0, 5i,_i,4 = 0, 1^2,-1, 



L\o = x'^ + y^ 



V\ 



o 



-2a 



iP + if 



25 



X 



26 



x(/3 + 7) 



n2 



25 



The set J is J = {(a,/3,7, (5) :5 < and 



-2a + 



iP+i? 



25 



> 0}. This set is feasible and 



(/3+7)' 



261 



the equations 5 = 0, —2a + 

request. 

Example 4. Let us have the system: 



0, with ^ < define the bifurcation variety upon 



X = —{y — a)x 

y = -ix^ 

where a and 7 are parameters taking values in R. The equilibrium points are (0, fc), A; G R. 
We shall try to find bifurcation varieties of the parameters for a specific equilibrium point 
(0, ^), ^ is a constant but otherwise arbitrary real number. We shall work with the following 
Lyapunov function candidate: 

L = Aix'^ + A2y'^ + A3X + A^y + A 

The Formal-Bif-Algorithm will give for L: 

L = A2{S2,-i,i + S2,i,ix + yf + (-2^2^2,1,1) (-^1,-1,2 + a^)(^2,-i,2 + ^2,1,22; + y) + 

+ (A - 2^^2,-1,1 + 2A25i,_i,2>S2,l,l) ("52,-1,3 + ^2,1,33: + y) + 
+ {Ai - ^2^2^,1 + 2^2^2,1,1^2,1,2) (^1,-1,4 + Xf + Ci(5i,_i,5 + x) + Rl 

For the quantity V , we get: 

V = (2Ai - 27yl2)(l^i,_i,i + xf{W2,-i,i + 1^2,1,1^: + y) + (A3 - 4AiW^i,_i,i+ 
+47A2m,_i,i)(m,_i,2 + a;) (1^2,-1,2 + ^^2,1,2^ + y) + {-2A^Wl_^^^ + 2-iA2Wl_^^^- 

-y43W"l,-l,2 + 4Aiiyi,_i,iiyi,_i,2 - 47A2iyi,_l,liyi,_l,2)(W^2,-l,3 + ^^2,133^ + ?/) + 
+ (-2^1^^2,1,1 + 27^2^^2,1,1) (W^l,-1,4 + xf + /ll(W^l,_l,5 + xf + /i2(W^l,-l,6 + X) + Rv 
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and thus the sets O and E are: 

O = {-2^2^2,1,1, A- 2^2^2,-1,1 +2A25i,_i,2^2,l,l, Ci, 2^-27^2, 

A3 - 4Ail^i,„i,i + 47^21^1,-1,1, . . . , /i2, ^^^2 + OAi + A5} 

£■ = {A2, Ai — ^25*2, 1,1 + 2yl.2S'2, 1,15*2, 1,2, hi, Rl, Rv} 

(The coefficients /ii,/i2 are large polynomial expressions of the parameters Wijk-) The 
nontrivial values of these parameters which make the members of the set O equal to zero 
are: 

= {Ai=-fA2, A4 = 2^2^2,-1,1, A, = -e^A2-9A,, ,5i,_i,4 = 0, 

^3 = 0, ^2,1,1=0, l^i,_i,5 = 0} 

and the quantities L and V become 

L = A2ix^ + ^2(^2,-1,1 + yf - A^ie + ^2,-1,1)' 

V = -2^ A2{a + 82,^1,1) x^ 

In order to ensure that the coefficients of the above polynomials are positive or equal to 
zero the next set must be feasible. 

J = {(Ai, A2, 7, e, 82,-1,1) ■■ A27 > 0, A2 > 0, -A2ie + ^2,-1,1)' > 0, -27yl2(a + ^2,-1,1) > 

The only non-trivial way to get that is: 

^2,-1,1 = -^, A2>0, 7 > 0, a + e<0 

and thus the bifurcation variety for the parameters a and 7 are {{a, 7) : a = —6, 7:7 = 0}. 
Example 5. The current example deals with a nonlinear system with three states: 

X = Lpy^ — xy + ax + 7y 

ij = x^ — tpxy + 5x + ipy 
z = -2z 
The equilibrium points are (0, 0, 0) and 

/-203 - 22av9 + aV^ ± 6^1764 + A55aip + 22a'ip'^ - a^v^^ 
V {7 + aip){13 + aif) ' 

-42 T v^T764T455a^T22^V^r^V ^ 







13(/9 + aip'^ y 

We shall work with the ffist one. We choose as Lyapunov function the classical one L 
x"^ + y"^ + z^ . The Formal Bif-Algorithm will give for V: 

V = 4(l^3,_i,i + l^3,i,iX + Ws,2,iy + Zf - 8W^3,2,l(W"2,-l,2 + W'2,1,2X + y\ 
■(1^3,-1,2 + W^,i,2X + Wz,2,2y + ^) - 8(1^3,1,1 " W"2,l,2V1^3,2,l) " (W"l,-1,3 + a:)" 
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■(W3,_i,3 + ^3,1,30; + 14/3,2,31/ + Z) - 8(m;3,_i,i - 1^1,-1,3^^3,1,1 - W^2,-1,2W^3,2,1 + 
+ W"l,-l,33W"2,l,2W"3,2,l)(W"3,-l,4 + 1^3,1,42: + W^3,2,4|/ + z) - 2{ip + 2Wl^^^ - 41^3,2,1^^3,2,2) " 

•(W"2,-i,5 + W2,i,5X + yy + ai(iyi,_i,6 + x) iW2,-i,e + W"2,i,6a^ + y)+ ^^2(1^2,-1,7 + 1^2,1,72: + y)+ 

+a3(iyi,_i,8 + xy + a4(W"i,-i,9 + x) + Rv 



(we do not write cti, o'2,cr3,cr4,i?y explicitly, due to their large size). The values O = { 

1^3,2,1 = 0, W^3,-l,l = 0, \ 
expressions for L and V: 



1^3,2,1 = 0, W^3,-i,i = 0, W2,-i,5 = 0, m,-i,8 = 0, 1^3,1,1 = 0, 11^2,1,5 = ^} will give the next 



r 2 I 2 I 2 

L = X + y + Z 

V\o=(-2a+—\x^-2^iy + —\ + Az'^ 

Obviously, for y? < 0, —2a + — > V^ is positive and hence the Lyapunov function is 
decreasing across the trajectories of the system and the origin is stable. Therefore, the 
relations </? = and —2a + — = 0, p < define the bifurcation variety upon request. 

5 Concluding Remarks 

The issue of this paper was the description of an algorithm which finds the so-called Hopf- 
Bifurcation varieties. These are relations among parameters, the violation of which changes 
the stability of an equilibrium point. This algorithm transforms polynomial expressions to 
sums of products of first degree polynomials with parametrical coefficients. By giving to 
these parameters proper values we ensure the positiveness of certain quantities and thus, 
we can investigate the stability behavior through Lyapunov theory. 
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